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We study the potential energy surface of the ozone molecule by means of Quantum Monte Carlo 
simulations based on the resonating valence bond concept. The trial wave function consists of an 
antisymmetrized geminal power arranged in a single-determinant that is multiplied by a Jastrow 
correlation factor. Whereas the determinantal part incorporates static correlation effects, the aug¬ 
mented real-space correlation factor accounts for the dynamics electron correlation. The accuracy 
of this approach is demonstrated by computing the potential energy surface for the ozone molecule 
in three vibrational states: symmetric, asymmetric and scissoring. We hnd that the employed wave 
function provides a detailed description of rather strongly-correlated multi-reference systems, which 
is in quantitative agreement with experiment. 

I. INTRODUCTION 


Beside its environmental and chemical importance, 
ozone plays a significant role in assessing the accu¬ 
racy of electronic structure methods due to its multi¬ 
reference character. The latter is due to static elec¬ 
tron correlation that arises in situations with degener¬ 
acy or near-degeneracy, as in transition metal chemistry 
and strongly-correlated systems in general. As a con¬ 
sequence, determining the electronic structure of ozone 
requires an accurate treatment of both, dynamic as well 
as static electron correlation effects. Hence, computing 
the total energy, equilibrium geometry and vibrational 
frequencies of the ozone molecule is a rigorous test for 
even the most accurate quantum-chem ical a pproaches, 
such as quantum Monte Carlo ( QMCl ^^, configura¬ 
tion interactiorP®, coupled-clustei®“ as we ll as multi- 
reference-based method^HMl, just to name a few. 

In this work, we study the potential energy surface 
(PES) of a single O 3 molecule using the resonating va¬ 
lence bond (RVB) theory of the chemical bond, which was 
first introduced by Pauling in the early days of quantum 
mechanic^ 37 . According to the RVB theory, a spin sin¬ 
glet can be formed between any two valence electrons be¬ 
longing to neighboring atoms. Within this approach, the 
ground state wave function of a molecule can lower the 
variational energy by imposing an auxiliary resonating 
degree of freedom between all valence bond configura¬ 
tions. However, the number of bonds increases exponen¬ 
tially with the number of atoms. Nevertheless, it was 
realized that a single-determinant wave function com¬ 
bined with a suitable real space correlation term, known 
as that Jastrow correlation fac toi^ , can be employed to 
represent a complex RVB stat d 39 * 4 ^ . To that extend, the 


variational energy is computed by means of conventional 
QMC techniques®®!, such as variationa l M onte Carlo 
(VMC) and diffusion Monte Carlo (DMCj^HS. The RVB 
wave function is expanded in terms of molecular orbitals 
(MOs) and optimized by a constrained energy minimiza¬ 
tion at constant number of I n this way, an ac¬ 

curate RVB wave function based on a single-determinant 
is obtained, which allows to study sys tems with sizable 
static and dynamic electron correlation ) 49 ! 50 ! 


II. RESONATING VALENCE BOND QUANTUM 
MONTE CARLO 

The RVB wave function, which is referred to as Jas¬ 
trow antisymmetrized geminal power (J-AGP), is given 
by a product of a Jastrow correlation factor J and a 
determinatal part d>jv, he. ^j-agp = Each 

geminal is described by a symmetric pairing-function 
= 4 > { r ji r i)- For an -/V-electron system, con¬ 
sisting of spin-up and AH spin-down electrons and 
assuming that AH < iW, we can define a pure spin state 
with a total spin S = |A"1 — |/2 and maximum spin 

projection S = , S* ot in terms of an antisymmetrized prod¬ 
uct of singlet pairs and 2 S unpaired orbitals ^(rj) of 

the remaining spin-up electrons. The resulting antisym¬ 
metrized geminal power (AGP) wave f unct ional, which 
can be arranged in a single-determinanlP 9 ®, reads as 

N l Af f 

^> N (r)=a n < K r i, T> i') n w 

i=1 j=AG + l 
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where R = |r|, • ■ ■ , rj^, rf, • • • , rj^j indicates the 

37V dimensional vector of all electron coordinates, 
while A is the antisymmetrization operator. In gen¬ 
eral, breaking number-symmetry leads to the Bardeen- 
Cooper-Schrieffer (BCS) wave function^, whereas break¬ 
ing spin-symmetry yields unrestricted Hartree-Fock 
determinantThis is to say that the RVB state corre¬ 
sponds to the Gutzwiller-projected BCS state, while the 
AGP is the particle-conserving version of the BCS wave 
function, which accounts for the static electron correla¬ 
tion. 

Dynamic correlation effects between the electrons is 
taken into account by the Jastrow correlation factor 


J(R) = exp 




. i<j 


( 2 ) 


where /(r^r^) is a two-electron coordinate function just 

as the previously introduced pairing-function 

and is chosen such as to satisfy the Kato cusp 

condition:^. 

Just like any other function of two coordinates, 
f( r i,rj) as well as c/)(ri,rj) can be expressed in terms 
of single-particle orbitals: 


The smallest number of orbitals, which is compatible 
with an unpolarized TV-electron system, is n = TVU In 
that case the rank of A is also minimal and identical 
to n and the antisymmetrization operator A singles out 
just one Slater determinant (SD). Due to the fact that 
the determinantal part of the resulting wave function is 
equivalent to the one of unrestricted Hartree-Fock (HF) 
theory, we will refer to it as JHF to emphasize the pres¬ 
ence of the Jastrow correlation factor. If the orbitals of 
the SD are determined by means of density functional 
theory (DFT), the eventual wave function is denoted as 
JDFT 57 . However, the great appeal of the AGP wave 
function is that if the rank of A is larger than TV7 , it cor¬ 
responds to a multireference wave function that is able 
to describe static electron correlation effects, although at 
the computational cost of a genuine single-determinant 
approach. We will name the resulting wave function as 
JAGP. As such, the RVB wave function can be thought 
of as the natural extension of the HF theory, where the 
Jastrow correlation factor accounts for the dynamic elec¬ 
tron correlation, while for n > 7W most of the static 
correlation is recovered. 


III. COMPUTATIONAL DETAILS 


M m 

4 > { r i: r j) = ( 3 ) 

a,b fJ.,v 

where tp(r) is an atomic orbital (AO), while to is the 
number of AOs per atom and M the number of atoms. 
The coefficients A“’j) are arranged in a quadratic m.M x 
rriM matrix denoted as A. Moreover, to conserve the 
total spin, A is also symmetric. Assuming that the AOs 
are non-orthogonal, i.e. = (<p| (fi) ^ I, diagonalizing 
A corresponds to solving the equation 

ASP = PA', (4) 

where the matrix P contains the eigenvectors of A and 
diagonal matrix A' the associated eigenvalues | | > 

I A!, | > ■ ■ ■ |A' mM | > 0. The number of nonzero eigen¬ 
values is equivalent to the rank of A. Using the identity 
P T SP = I that is identical to P T = (SP ) -1 and multi¬ 
plying it to the right of Eq. Q, A = PA' P T is obtained. 
From this it follows that Eq. ([3]) is equivalent to 

n—2S 

<f>i r iAj)= X 'kM r i)M r j) ( 5 ) 

k =1 

where n is the number of MOs (j>k{fi) that can be written 
as a linear combination of AOs: 

fain) = ^2c k m(ri), ( 6 ) 

i 

where Cki are the AO coefficients that independently on 
the particular basis set used are implicitly dependent on 

Vi- 


All of our QMC calculations were performed using the 
TurboRVB program package of Sorella and coworkerd^. 
We have employed the two different trial wave functions, 
the JDFT and the JAGP wave function. In both cases 
the single-determinant is made of MOs as obtained by 
DFT within local-density approximation (LDA), as de¬ 
scribed in Ref. EB In case of the JAGP wave func¬ 
tion all parameters were variationally optimized at the 
VMC level of theory using the stochastic reconfigura¬ 
tion algorithm®®!!^ while for the JDFT wave function 
the parameters of the Jastrow factor were held fixed. 
Scalar-relativistic energy consistent pseudopotentials— : 
were used to describe two core electrons of the oxy¬ 
gen atoms. For the single-determinant we employed an 
atomic basis set consisting of contracted 10s, 8p, 6 d, 4/ 
Gaussian-type orbitals, while for the Jastrow correla¬ 
tion factor uncontracted 2s, 2 p, 1 d Gaussian-type orbitals 
were used. The variational energies were computed us¬ 
ing two different QMC methods, VMC and the even 
more accurate lattice-regularised diffusion Monte Carlo 
(LRDMC) projection teclmiqutP^, which eventually re¬ 
sults in four trial wave functions we have investigated 
here: JDFT-VMC, JDFT-LRDMC, JAGP-VMC and 
JAGP-LRDMC. 


IV. RESULTS AND DISCUSSION 

The PES of the 0 3 molecule exhibits two types of min¬ 
ima, one of C 2 V symmetry (open minimum) with an apex 
angle of 116.75° and another one of D 3 h symmetry (ring 
minimum), where the oxygen atoms form an equilateral 
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FIG. 1. (colour online) The PES of the Os(X) molecule in a 
symmetrical vibration state as obtained by the JDFT-VMC, 
JDFT-LRDMC, JAGP-VMC and JAGP-LRDMC methods. 
The energy normalization is such that zero energy corresponds 
to the ground state energy of 02( 3 Xjj") + 0( 3 P). The experi¬ 
mental dissociation energy is illustrated by the dotted line. 


FIG. 2. (colour online) The PES of the Oa(X) molecule in an 
asymmetrical vibration state as obtained by the JDFT-VMC, 
JAGP-VMC, JDFT-LRDMC and JAGP-LRDMC methods. 
The energy normalization is such that zero energy corresponds 
to the ground state energy of 02( 3 £j ) + 0( 3 P). The experi¬ 
mental dissociation energy is illustrated by the dotted line. 


triangle. Due to the fact that the open structure is en¬ 
ergetically significantly lower than the closed ring case, 
we will confine ourselves to the O 3 molecule with C 2 „ 
symmetry only. Symmetry considerations suggest the 
existence of three open minima, which are separated by 
energetically high barriers and thus can be considered as 
independent from each othei®. As a consequence, the 
potential energy E(d±,d 2 ,a) is calculated as a function 
of just three variables, where d\ and c ?2 are the two 0 c —0 
bond distances, while a is the O — O c — O bond angle 
and O c the central atom. The experimental values for d\, 
d 2 , a, and the dissociation energy are 2.4052ao, 2.4052oo, 
116.75°, and 26.105 kcal/mol , respectivel}^. 

Figure [l] illustrates the PES of the Os(X) molecule 
in a symmetrical vibration state, where the value of 
a is fixed at 116.8° and the d\ = d 2 - The com¬ 
pressed 0 3 (X) molecule is dominated by Coulomb re¬ 
pulsion and the results of all trial wave functions in¬ 
vestigated here rather similar. However, whenever the 
0 3 (X) molecule is stretched, static electron correlation 
becomes more important and the JAGP wave function 
superior due to the strong multi-reference character of 
the stretched O — O bond. The dissociation energies at 
the experimental equilibrium geometry, using the JDFT- 
VMC, JDFT-LRDMC, JAGP-VMC, as well as JAGP- 
LRDMC methods are -16.60(5), -19.74(4), -25.25(6) and 
-26.14(5) kcal/mol, respectively. In fact, even in at the 
variational level, the JAGP wave function outperforms 
accurate LRDMC calculations based on the JDFT trial 
wave function and recovers ~ 97% of the experimental 
dissociation energy, while at the JAGP-LRDMC level the 
result is correct up to 0.1%. 

The PES of the Os(X) molecule in an asymmetrical 
stretching state is shown in Figure [2j where the value 
of a is again fixed at 116.8°. Defining Sd as deviation 



Angle (degree) 

FIG. 3. (colour online) The PES of the Oa(X) molecule in 
a scissoring vibration state as obtained by the JDFT-VMC, 
JDFT-LRDMC, JAGP-VMC and JAGP-LRDMC methods. 
The energy normalization is such that zero energy corresponds 
to the ground state energy of 02( 3 S ! j') + 0( 3 P). The experi¬ 
mental dissociation energy is illustrated by the dotted line. 


from the equilibrium bond length, each point on the PES 
corresponds to d\ ± Sd and d 2 T Sd. Given the symme¬ 
try, E(di + Sd, c ?2 — Sd,a) = E{d\ — Sd,d 2 + 6 d,a), as 
can been seen in the axial symmetry around the equilib¬ 
rium bond length in Figure [2] The corresponding dis¬ 
sociation energies for the JDFT-VMC, JDFT-LRDMC, 
JAGP-VMC, as well as JAGP-LRDMC techniques are 
-16.58(6), -19.71(5), -24.99(7) and -25.64(5) kcal/mol, 
respectively. Interestingly, at equilibrium the JAGP is 
significantly superior to the JDFT trail wave function 
regardless of the particular QMC method, whereas for 
large asymmetries the LRDMC technique is essential. 
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Figure [i] exhibits the PES of the O 3 (X) molecule 
in a scissoring state, where d\ and d 2 are both fixed 
at 2.4052ao and, thus, a the only variable. At the 
equilibrium bond angle, the dissociation energies as ob¬ 
tained by the JDFT-VMC, JDFT-LRDMC, JAGP-VMC 
and JAGP-LRDMC methods are -17.42(7), -20.55(5), - 
25.61(6), and -25.97(5) kcal/mol, respectively. Again, 
the JAGP systematically improves the description of 
the JDFT wave function. In fact, even at the JAGP- 
VMC method exceeds JDFT-based LRDMC calcula¬ 
tions. Moreover, the small difference between the JAGP- 
VMC and JAGP-LRDMC schemes clearly demonstrates 
accuracy of the JAGP approach. 


all cases we have considered here. In particular, ex¬ 
cept for the most asymmetric configurations, the JAGP- 
VMC method was able to outperform even more accurate 
LRDMC calculations based on the JDFT trail wave func¬ 
tion. Together with the fact that the differences between 
the JAGP-VMC and JAGP-LRDMC results were rather 
small, this is a clear manifestation of the superiority of 
the JAGP wave function. In any case, the significant 
deviation between the JDFT and JAGP trail wave func¬ 
tions highlights that for a multi-reference system, such 
as the 0 3 (X) molecule, an accurate treatment of static 
electron correlation effects is essential. 


V. CONCLUSIONS 

In summary, we find that the JAGP systematically im¬ 
proves upon the conventional JDFT wave function in 
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